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We have carried out a detailed and systematic study of the correlation energies of inert gas atoms 
Ne, Ar, Kr and Xe using relativistic many-body perturbation theory and relativistic coupled-cluster 
theory. In the relativistic coupled-cluster calculations, we implement perturbative triples and include 
these in the correlation energy calculations. We then calculate the dipole polarizability of the ground 
states using perturbed coupled-cluster theory. 
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I. INTRODUCTION 

High precision atomic experiments are at the core of 
several investigations into fundamental physics and high 
end technology developments. Selected examples are 
search for electric dipole moment (EDM) jT] and obser- 
vation of parity nonconservation [2]. These endeavours, 
in general, require precision atomic theory calculations 
to analyse the results and understand systematics. The 
challenging part of precision atomic structure and prop- 
erties calculations is obtaining accurate wave functions. 
In the case of high Z atoms, the need to incorporate 
relativity adds to the difficulty. A systematic study of 
the correlation energy is one of the possible methods to 
test the accuracy of the atomic wave function. In this 
paper, we report the results of correlation energy calcu- 
lations of inert gas atoms Ne, Ar, Kr and Xe. For this 
we employ many body perturbation theory (MBPT) and 
calculate the second order correlation energy. A com- 
parative study reveals the changing nature of electron 
correlations in the group. Our interest in particular is 
Xe, which is a candidate for EDM experiments [3] and 
theoretical calculations |4 . 

For completeness, in the presentation of the paper, we 
give an overview of MBPT. It is a powerful theory and 
forms the basis of other more sophisticated and elaborate 
many-body methods. However, one drawback of MBPT 
is the complexity of expressions at higher orders. This 
renders the theory inappropriate to incorporate strong 
correlation effects in heavy atoms. Yet, at lower orders 
its simplicity makes it an ideal choice to test and optimize 
basis sets. We use this insight to generate basis sets for 
coupled-cluster calculations. 

The coupled-cluster theory, first developed in nuclear 
many body physics [5j [6], is considered the most accu- 
rate many body theory. In recent times, it has been used 
with great success in nuclear [7], atomic [HI IS], molec- 
ular [TU] and condensed matter [IT] calculations. It is 
equivalent to incorporating electron correlation effects to 
all orders in perturbation. The theory has been used 
in performing high precision calculations to study the 
atomic structure and properties. These include atomic 
electric dipole moments [H [12], parity nonconservation 
[IB] , hyperfine structure constants [HI HI] and electro- 



magnetic transition properties |15[ 116] . In the present 
work we use the relativistic coupled-cluster singles and 
doubles (CCSD) approximation to calculate correlation 
energy and dipole polarizability of inert gas atoms Ne, 
Ar, Kr and Xe. In the dipole polarizability calculations, 
the dipole interaction Hamiltonian is introduced as a per- 
turbation. A modified theory, recently developed |17j . 
incorporates the perturbation within the coupled-cluster 
theory. This theory has the advantage of subsuming cor- 
relation effects more accurately. The results provide a 
stringent test on the quality of the wave functions as 
the dipole polarizability of inert gas atoms are known 
to high accuracy [18] ■ Based on the CCSD method, we 
also estimate the third order correlation energy. Fur- 
ther more, perturbative triples are incorporated in the 
coupled-cluster calculations. 

In the paper we give a brief description of MBPT in 
Section. II and discuss the method to calculate electron 
correlation energy to the second and third order in resid- 
ual Coulomb interaction. The coupled-cluster theory is 
described in Section. Ill, where we also discuss linearized 
coupled-cluster theory and correlation energy calculation 
using coupled-cluster theory. Then the inclusion of ap- 
proximate triples to the correlation energy is explained 
and illustrated. Section. IV is a condensed description of 
the perturbed coupled-cluster theory and provide details 
of how to incorporate the effects of an additional pertur- 
bation to the residual Coulomb interaction in atomic sys- 
tems. Results are presented and discussed in Section. V. 
In the paper, all the calculations and mathematical ex- 
pressions are in atomic units (e = h = m e = 1). 



II. CORRELATION ENERGY FROM MBPT 

In this section, to illustrate the stages of our calcula- 
tions and compare with coupled-cluster theory, we pro- 
vide a brief description of many-body perturbation the- 
ory. Detailed and complete exposition of the method, in 
the context of atomic many-body theory, can be found 
in ref [15] , 

The Dirac-Coulomb Hamiltonian H BC is an appro- 
priate choice to incorporate relativistic effects in atoms. 
This is particularly true for heavy atoms, where the rel- 
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ativistic effects are large for the inner core electrons due 
to the high nuclear charge. As the name indicates, H DC 
is fully relativistic for the one-body terms only. For an 
N electron atom 

N 1 
H DC = ^ [c«, • Pi + ((3- l)c 2 - V N (n)] + ^ — , 

i ■ , ■ r ij 

2—1 1<J J 

where on and /3 are the Dirac matrices. For the nuclear 
potential Vpf(r), we consider the nucleus as finite size and 
the volume effects are taken into account by modeling 
the nuclear charge distribution as a Fermi two-parameter 
distribution. Then the nuclear density is 

Pnuc(r) = (2) 

here, a — i41n3. The parameter c is the half-charge 
radius, that is p n uc(c) = po/2 and t is the skin thick- 
ness. The eigen states of £f DC are the correlated 
many-particle states with eigenvalues Ei. The eigenvalue 
equation is 

if DC |^)=^|^). (3) 

It is however impossible to solve this equation exactly 
due to the relative coordinates in the electron-electron 
Coulomb interaction. Many-body perturbation theory 
is one approach which, starting from a mean field ap- 
proximation, incorporates the electron correlation effects 
systematically. 

The starting point of perturbative scheme in many- 
body theory is to split the Hamiltonian as 

H DG = H + V, (4) 

where H = J2d ca i-Pi + (A ~ l )° 2 ~ v n^ + w(n)], is 
the unperturbed or zeroth order Hamiltonian. It is the 
exactly solvable part of the total Hamiltonian and cor- 
respond to independent particle model. In this model, 
each electron is assumed to move independently of the 
others in an average field arising from the nucleus and 
other electrons. The average field of the other electrons 
is the Dirac- Fock central potential u(r^). The remain- 
ing part of the electron-electron Coulomb interaction 
V = J2i<j T~ ~ J2i u ( r i)> I s the residual Coulomb in- 
teraction. The purpose of any atomic many-body theory 
is to account for this part as accurately as possible. The 
Hamiltonian Hq satisfies the eigenvalue equation 

H \$i) = E$\$i), (5) 

where is a many-particle state and Ef is the eigen- 
value. The cigcnstatcs are generally Slater determi- 
nants, antysymmetrised direct product of single particle 
states and Ef is the sum of the single particle energies. 
The difference between the exact and mean field energy, 
AEi = Ei — E®, is the correlation energy of the i Th state. 



At the single particle level, the relativistic spin orbitals 
are of the form 

ib (r)--( P ™( r )x K m(r/r) \ 

where P nti (r) and Q nK (r) are the large and small com- 
ponent radial wave functions, k is the relativistic total 
angular momentum quantum number and X K m( r / r ) are 
the spin or spherical harmonics. One representation of 
the radial components is to define these as linear com- 
bination of Gaussian like functions and are referred to 
as Gaussian type orbitals (GTOs). Then, the large and 
small components [501 HI] are 

p 

Q nK {r) =J2C s Kp g s Kp (r). (7) 
p 

The index p varies over the number of the basis functions. 
For large component we choose 

SkpW = ^r"-e-^ r \ (8) 

here n K is an integer. Similarly, the small component are 
derived from the large components using kinetic balance 
condition. The exponents in the above expression follow 
the general relation 

ap^aopP- 1 . (9) 

The parameters ao and (3 are optimized for an atom to 
provide good description of the atomic properties. In 
our case the optimization is to reproduce the numerical 
result of the total and orbital energies. Besides GTO, 
B-splines is another class of basis functions widely used 
in relativistic atomic many body calculations [22]. A 
description of B-splines with details of implementation 
and examples are given in ref |23j . 

The next step in perturbative calculations is to divide 
the entire Hilbert space of H into two manifolds: model 
and complementary spaces P and Q respectively. The 
model space has, in single reference calculation, the eigen 
state of H which is a good approximation of the 
exact eigenstate |^,) to be calculated. The other eigen- 
states constitute the complementary space. The corre- 
sponding projection operators are defined as 

P=|$i>($i| and Q= ^ |*i><*j|. (10) 

The operator P projects out the component of the exact 
eigenstate which lies in the model space, P\^i) = l^i) 
and Q projects out the component in the orthogonal 
space and P+Q = 1. In the present paper, we restrict to 
calculating the ground state \^o) of the closed shell inert 
gas atoms. From here on, for a consistent description, 
the model space consist of |$o)- 
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The most crucial part of perturbation theory is to de- 
fine a wave operator f2 which operates on |<f>o) and trans- 
form it to |W ) as 



(11) 



Then, with the intermediate normalization approxima- 
tion (Wq|$o) = 1; the wave operator is evaluated in or- 
ders of the perturbation as = Y,Zo ftW with ni ° } = L 
It is possible to evaluate fiW iteratively or recursively 
from the Bloch equation 



[n, h }p = Qvnp - xPVQP, 



(12) 



where x — 1S the correlation operator. In the 

second quantization notations, the wave operator and 
perturbation can be expressed in terms of particle ex- 
citations. Then, the effect of correlation is incorporated 
as linear combination of excited states. For simplifica- 
tion, in the normal form the perturbation is separated as 
[Hi] V = Vo + V\ + V 2 . These zero-, one- and two-body 
operators and are defined as 



Vi 

v 2 



core 1 core 

= 5>| - «|o) + - J2« ab Kb\ ab ) H^lofc)). 



ab 



= 53{ ffl i0 3 -}(i|«b'>, 

ij 

= ~^2{ a l a ] a i a k}(ij\ri2\ kl )- 



(13) 



ijkl 



The operators cr (a) create (annihilate) electrons in vir- 
tual (p,q,r,s, .... etc) and core (a,b,c,d, .... etc) shells. 
The indexes k, I etc are general representations of or- 
bitals, it could either be core or virtual. The operator 
Vo acting on |<j> ) leaves it unchanged, while V\ and V 2 
produce single and double excitations. From these def- 
inition, the first-order wave operator can be separated 
as 



n« =nf» + n 



(i) 



(14) 



Here, and Q% are one- and two-body components 
of the first order wave operator. The corresponding al- 
gebraic expressions are 



(i) 



ap ta e P 

o(i) ! v t t (j>q\v\ab) 

$2 2 = 9 2^ a p a q a b a a- 



abpq 



We get singly (doubly) excited states ( 

ill (Hi ) operates on the reference state 
plcxity of the expressions increases with 
turbation and is hard to manage. One 
in many-body perturbation theory is the 
evaluation of the perturbation expansion, 
dious algebraic evaluations are reduced to 



(15) 
(16) 



\<I>ab pq )) when 

$o . The corn- 
order of per- 
powerful tool 
diagrammatic 
Then, the te- 
a sequence of 



diagrams and equivalent algebraic expressions are derived 
with simple rules. Even with this approach, it is compu- 
tationally not practical to go beyond fourth order. 



n-2 



FIG. 1: The diagrammatic representations of the one- and 
two-body wave operator. Lines with downward (upward) ar- 
rows represent core (virtual) single particle states. 



A. Second-Order Correlation Energy 

The ground state correlation energy AEq, as described 
earlier, is the difference between the exact energy and 
the mean field energy. It is the sum total of the energy 
corrections from all orders in perturbation. At the n th 



order, the energy correction E, 



(n) 



and AE = £, 
energy is 



E 



(n) 



(*o|^ (n_1) |*o> 



Then the second order correlation 



Eill = (<f> \(v 1 + v 2 )(ny + n^)\<s> 



(17) 



When Dirac-Fock orbitals are used, the diagonal matrix 
elements of V\ are the orbital energies and off diagonal 
matrix elements arc zero. For this reason, it does not 
contribute to the second-order energy. Then, the second 
order correlation energy is 



E£l = ($ n \VMi 1) \* 



(18) 



There are two diagrams arising from the above expression 
and these are as shown in Fig(2j 





FIG. 2: MBPT diagrams arising from the second-order cor- 
relation energy. 



In terms of algebraic expressions 



corr / j 
abpq 



(ab\ V 2 \pq) (pq\V 2 \ab) (ba\V 2 \pq) (pq\V 2 \ab) 



— e g e a + £b — e P — e 9 

(19) 

In the above expression, the first and second terms on the 
right hand side are the direct and exchange. Though the 
expression is fairly straight forward to derive, we have 
given explicitly for easy reference while analysing the re- 
sults. 



B. Third-Order Correlation Energy 

The diagrammatic representation of the second order 
wave operator consists of single, double, triple and 
quadruple excitations. The singles are non-zero starting 
from the second order when Dirac-Fock orbitals are used. 
And, the triples and quadruples begin to contribute from 
this order. The triples consist of connected diagrams, 
whereas all the quadruples are disconnected. The third 
order correlation energy is 

Egl = ^ \(v 1+ v 2 )(n^+ni 2) +n^+n{ 2) )\^). (20) 

The triple and quadruple excitations do not contribute as 
V at the most can contract with double excitations. For 
the same reason mentioned earlier, in second order, V\ 
also does not contribute. Then the third order correlation 
energy is simplified to 



($o|y 2 ^ 2) ) 



l*o) 



(21) 



This is similar in form to the second order correlation 
energy. In general, the n th order correlation energy has 
non-zero contribution from the term V^fi^™ only. It 
must be mentioned that the connected triples begin to 
contribute from the fourth order energy This is utilized 
in perturbative inclusion of triples, in later sections of 
the paper, while discussing coupled-cluster calculations. 



n 



(0) 



t: 
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FIG. 3: Diagrammatic representation of unperturbed single 
and double cluster operators. 



cluster operators in the second quantized notations are 
7\ (0) = £i£4a a , (24) 



r (0) _ I V t pq <Ja 
2 — 2! ^— J p 



\a b a a . 



(25) 



a,b,p,q 



Here, t p and t P l are the single and double cluster am- 
plitudes respectively. Subtracting (^ol-H" |<I>o) from both 
sides of Eq.d3|) and using the normal form of an operator, 
N = O - (3» |0|$o), we get 



^n|*o) 



0/ 



(26) 



where AE = E — (<5> |-ff l^o), as defined earlier, is the 
correlation energy. Operating with e~ T and projecting 
the above equation on excited states we get the cluster 
amplitude equations 



(*SI^_n|*o) = 0, 

M#n|<I>o)=0, 



(27) 
(28) 



III. COUPLED-CLUSTER THEORY 

The coupled-cluster theory is a non-pcrturbative 
many-body theory and considered as one of the best. A 
recent review [53] provides an excellent overview of recent 
developments and different variations. In the context of 
diagrammatic analysis of many-body perturbation the- 
ory, coupled-cluster theory is equivalent to a selective 
evaluation of the connected diagrams to all orders. Then 
casting the disconnected but linked diagrams as products 
of connected diagrams. In coupled-cluster theory, for a 
closed-shell atom, the exact ground state is 



l*o> 



(22) 



where = e _T< ' H-^e T< ' '' is the similarity transformed 
or dressed Hamiltonian. Following Wick's theorem and 
structure of An, in general 



H N = H N + {&>} + I{F N T(°)T(°)} 



i 1 

Here A . . . B denote contraction between two operators 
A and B. The single and double cluster amplitudes are 
solutions of Eq.(27l and (28 1 respectively. These are set 



of coupled nonlinear equations and iterative methods are 
the ideal choice to solve these equations. 



where is the cluster operator. The superscript is 
a tag to identify cluster operators arising from different 
perturbations. For the case of N electron atoms, the 
cluster operator is 



N 



T (0) 



i=l 



(0) 



(23) 



In closed shell atoms, the single and doubles provide a 
good approximation of the exact ground state. Then, 
the cluster operator — + and is referred to 
as the coupled-cluster single and doubles (CCSD). The 



A. Linearized Coupled-Cluster 

The nonlincarity in the cluster amplitude equation 
arises from the two and higher contractions in the dressed 
Hamiltonian. An approximation often used as a starting 
point of coupled-cluster calculations is to retain only the 
first two terms in H n, then 



H N = H N + {H N T^}. 



(30) 
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FIG. 4: Diagrams which contribute to the singles, unper- 



turbed cluster operator (T x 
cluster. 



(OK 



in the linearised coupled- 



B. Correlation Energy and approximate triples 



From Eq. ( 26 ) the ground state correlation energy, in 
coupled-cluster theory, is the ground state expectation 
value of iJjy- That is 



AE= ($ |ff N |$o) 



(35) 



The diagrams arising from the above expression arc 
shown in the Fig|6] The dominant contributions are from 
the diagrams (a) and (b), which is natural as the doubles 
cluster amplitudes are larger in value than the singles. 
Diagram (e) does not contribute to the correlation en- 
ergy when Dirac-Fock orbitals are used. 









o 











(b) 



(c) 



(f) 



fe) 



(li) 





FIG. 6: Coupled-cluster correlation energy diagrams. The 
diagram (e) is equal to zero when Dirac-Fock orbitals are used. 



FIG. 5: Diagrams which contribute to the doubles, unper- 
turbed cluster operator (Tj ), m the linearised coupled- 
cluster. 



The cluster equations are then a pair of linear equations 



<*S|J3- N + {JT N r< >}|$o) = 0, 
«|# N -HiS (0) }|<i>o) = o. 



(31) 
(32) 



In the CCSD approximation T<W = T 1 (0) + T 2 W , these 
equations are then 

(^{iSf} + {iSr 2 (0) }|$ > = -(*g|ff N |*o) 

«|{S 0) } + {S 2 (0) }|<I>o) = -(Kt\H^o). 

(33) 

These are the linearized coupled-cluster equations of sin- 
gle and double cluster amplitudes. This can be combined 
as the matrix equation 



,(0) 



H21 H22 



H20 



(34) 



where H n = (*S|fT N |$g), H 12 = ($£|.ff N |$£) and so 
on. The equations are set of coupled linear equations 
and solved using standard or specialized linear algebra 
solvers. In the literature several authors refer to lin- 
earized coupled-cluster as all-order method. A descrip- 
tion of the all-order method and applications are given in 
ref |25j . In a recent work, the authors report the combi- 
nation of all-order method and configuration interaction 



To go beyond the CCSD approximation, we incorpo- 
rate selected correlation energy diagrams arising from ap- 
proximate triples. The approximate triples are perturba- 
tive contraction of V 2 with the cluster amplitudes 
[2"T1 [2"B] . Example diagrams of the approximate triples 
and correlation energies are shown in Fig|7] There are 
two categories of triples, first is V2 contracted with 
through a hole line, and second contraction through a 
particle line (Fig{7^,). To calculate the correlation energy 
from the triples, these are contracted perturbatively with 
V2 and reduced to a double excitation diagram. Then 
the correlation energy is obtained after another contrac- 
tion with V 2 . These two contractions generate several 
diagrams. The triples correlation energy diagrams are 
separated into three categories based on the number of 
internal lines. These are: two particle and two hole in- 
ternal lines (2p-2h), three particle one hole internal lines 
(3p-lh), and one particle three hole internal lines (lp- 
3h). In the present calculations eight diagrams from the 
first category and two each from other remaining two 
categories are considered. 



IV. PERTURBED COUPLED-CLUSTER 

The atomic properties of interest are, in general, asso- 
ciated with additional interactions. The interaction are 
either internal like hyperfrne interaction or external like 
static electric field. These are treated as perturbations 
which modify the wave function and energy of the atom. 
This section briefly describes a method to incorporate 
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The dressed Hamiltonian H^ is same as in Eq. ( 29 1 . Like 



in linearized coupled-cluster, these form a set of linear 
algebraic equations. 



FIG. 7: Diagrams of approximate triples calculated pertur- 
batively:(a) approximate triples cluster operator and (b) cor- 
relation energy arising from approximate triples. 



n(l) 



t: 



(i) 



FIG. 8: Diagrams of single and double perturbed cluster op- 
erators. 



an additional perturbation within the frame work of rel- 
ativistic coupled-cluster. The scheme is referred to as 
perturbed coupled-cluster theory. It has been tried and 
tested in precision atomic properties and structure calcu- 
lations. In the presence of a perturbation Hi, the eigen 
value equation is 



DC 



XHi) |* > =E\* ). 



(36) 



Here |^o) is the perturbed wave function, E is the corre- 
sponding eigenvalue and A is the perturbation parameter. 
The perturbed wave function is the sum of the unper- 
turbed wave function and a correction |^ ) arising from 
Hi. That is 



|*o) = |*o) + A|* c 



(37) 



Following the earlier description of coupled-cluster wave 
function, the perturbed wave function is 



l*o> 



l*o>. 



(38) 



The cluster operators T^°\ as defined earlier, incorporate 
the effects of residual Coulomb interaction. For clarity 
these are referred as unperturbed cluster operator. The 
cluster operators arise from H\ and are referred to 
as the perturbed cluster operators. It acts on the refer- 
ence state |$o) t° generate the correction. Consider the 
perturbation expansion to first order in A, we get 



|$o) = e T<0, (l 



ATW)|*o) 



(39) 



then 



To derive the cluster equations use this in Eq.(36l 
operate with e _T(0> and project on excited states. We get 
the equations for singles and doubles perturbed cluster 
amplitudes 



(^ a \{H N T^}\<i> Q ) 

(*^|{£r (1) }|*o) 



(Kl\Hi\^o) 



(40) 
(41) 



A. Approximate triples 

Like in T(°\ a perturbed triple cluster Fig(9|a) is a 

perturbative contraction between V2 and T^ 1 . As in the 
case of unperturbed approximate triples discussed ear- 
lier, there are two types of diagrams in the present case 
as well. One arises from particle line contraction and the 
other from hole line contraction between V2 and T^ 1 ' di- 
agrams. In this work we implement approximate triples 
while calculating properties. In particular, to calculate 
dipole polarizability and an example diagram is shown in 
Fig. [9jb). The algebraic expression of this diagram is 

\- {ab\4 0)i \pq){c\d\ S ){qs\V 2 \rc){pr\TP\ab) 



a,b,c,p,q,r,s 



(42) 

here, d is the dipole operator. In total there are twen- 
tyfour properties diagrams arising from the perturbative 
triples and we include all of these diagrams in the calcu- 
lations. 



B. Dipole Polarizability 

When an atom is placed in an external electric field £ , 
the charge distribution of electron cloud is distorted and 
an electric dipole moment D m d is induced. The dipole 
polarizability of the atom a is then the ratio of the in- 
duced dipole moment to the applied electric field, that 
is 



A nd = a£. 



(43) 



By definition, the dipole polarizability of the ground state 
is 



-2E 



\(*o\D\*i) 

E - El 



(44) 



where are intermediate atomic states. These are op- 
posite in parity to the ground state |^o)- The expression 
of a can be rewritten as 



a= -2(9 \D\9 ) 



(45) 



Here |¥„) = E/(l*i)(*i|-D|*o))/(^o - Ej), which fol- 
lows from the first order time independent perturbation 
theory. The perturbation Hamiltonian is Hi — —D.£ 
and external field £ is the perturbation parameter. One 
short coming of calculating a from Eq.(44l is, for practi- 
cal reasons, the summation over / is limited to the most 
dominant intermediate states. However, the summation 
is avoided altogether when the perturbed coupled-cluster 
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FIG. 9: Diagrams of approximate triples calculated per- 
turbatively:(a) (a) Representation of approximate perturbed 
triples, (b) Contribution of approximate perturbed triples to 
the dipole polarizability. 



wave functions are used in the calculations. From Eq. ( 38 ) 
the perturbed wave function 



I*, 



(46) 



In a more compact form, the dipole polarizability in 
terms of the perturbed coupled-cluster wave function is 



o 



(*o|£>|*o) 



(47) 



After simplification, using the perturbed wave function 



in Eq. (38 1, we get 



a= (* |D|* ) + (* |r>|*c 



(48) 



The correction I^q), as described earlier, is opposite in 



parity to |$ )- 



Hence the matrix elements (Wol-D^o) 
and (^qIUI^q) are zero. As D is hermitian, the two 
terms on the left hand side are identical and the above 
expression is same as Eq.(45l. Considering the leading 
terms 



(49) 



Here, the operator = e T<0) 



De Ti ' is the unitary 



transformed electric dipole operator. It is explicitly evi- 
dent that the dipole polarizability, in terms of perturbed 
cluster operator, does not have a sum over states. In this 
scheme, contributions from all intermediate states within 
the chosen configuration space are included. For preci- 
sion calculations, this is a very important advantage. 



V. RESULTS 

In order to get reliable results, in atomic structure and 
properties calculations, one prerequisite is good quality 
orbital basis set. In all calculations described in the pa- 
per, we employ GTOs as orbital functions. In particular, 
we use even tempered basis in which the parameters cto 
and (3, in Eq. ([9|, are different for each symmetries. We 
use the basis parameters of Tatewaki and Watanabe [25] 

p(0) 



energies of core orbitals) in agreement with the numeri- 
cal results. The numerical results are obtained from the 
GRASP92 30j code. In order to obtain converg ed E&ir, 
we consider orbital basis set consisting of all the core or- 
bitals and virtual orbitals up to 10,000-11,000 in single 
particle energies. 

The working equations of coupled-cluster theory are 
coupled nonlinear equations. Solving these equations is 
a computational challenge. The number of unknowns, 
cluster amplitudes, are in the order of millions. In addi- 
tion, implementing fast and efficient algorithms demand 
huge memory to tabulate and store two-electron inte- 
grals. This is essential as the two-electron integrals are 
needed repeatedly and are compute intensive. For the 
larger basis sets in the present work, the number of two- 
electron integrals is more than 2 x 10 s . In order to utilize 
memory efficiently, we have developed a scheme which 
parallelize the tabulation and storage of two-electron in- 
tegrals. To improve performance further, we also tabu- 
late and store 6j symbols. This is desirable as the angu- 
lar part of the perturbed cluster diagrams involve large 
number of 6j symbols. Quantitatively, we observe a per- 
formance gain of 30% or more with the 6j symbols tab- 
ulation. We shall describe and discuss the various com- 
putational schemes developed in a future publication. 

The unperturbed and perturbed cluster amplitude 
equations are solved iteratively using Jacobi method. We 
chose the method as it is relatively simple to parallelize. 
One drawback of the method is, it converges slowly. To 
obtain faster convergence, we employ direct inversion in 
iterated subspace (DIIS) [31] convergence acceleration. 



A. Second-Order Correlation Energy 

The SCF energy E^q, second-order correlation energy 

(2) 

Ecorr and the total energy E from our calculations are 
listed in Table |TJ For comparison the results of previ- 
ous calculations arc also listed. It is evident that our 
second-order correlation energy and total energy, sum 
of the SCF and second order correlation energy, are in 
agreement with the results of Ishikawa et al [32] for all 
the atoms studied. Among the previous works, we select 
the results of Ishikawa et al [32] for detailed comparison 
as their calculations are relativistic. The other results 
listed in the table are from non-relativistic calculations. 
For all the atoms, our SCF energy Eq C are lower and 
there are no discernible trends, as a function of nuclear 
charge, in the difference. Interestingly, except for Xe, 
our second order correlation energies are higher. This 
compensates the lower E^q and subsequently, the total 
energies E of the two calculations are in excellent agree- 
ment. The lack of trend indicates the choice of the basis 
set could be the source of the observed differences of E^X 

and Ecorr between the two calculations. The Table II 



as starting values and optimized further to obtain Efcfc lists the cumulative contributions from various symme- 
(ground state Dirac-Fock energy) and e a (single particle tries to Ecorr- Among the previous works, Lindgren and 
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TABLE I: The SCF E^ c , the second-order correlation 
listed are in atomic units (Hartrees). 



and the total energies E of Ne, Ar, Kr and Xe. All the values 



Atom 


Z 


Atomic mass 




This work 






Other work 










^DC 


E {2) 

-^corr 


E 


^DC 




E 


Ne 


10 


20.18 


-128.6932 


-0.3830 


-129.0762 


-128.6919 


-0.3834" 
-0.3836'' 
-0.3822 c 
-0.3697 d 
-0.3804 e 


-129.0753 


Ar 


18 


39.95 


-528.6882 


-0.6938 


-529.3820 


-528.6838 


-0.6981" 
-0.6822 6 
-0.685 / 
-0.790 9 


-529.3819 


Kr 


36 


83.80 


-2788.8659 


-1.8426 


-2790.7085 


-2788.8615 


-1.8468" 


-2790.7083 


Xe 


54 


131.29 


-7446.8887 


-2.9767 


-7449.8654 


-7446.8880 


-2.9587" 


-7449.8467 



"Referenced 
b Reference |33 
c Reference [3^ 
d Reference |36 
e Referencel37 
■^Reference [38 
s Reference|34 



1.5 r 



1.0 ■ 



a 0.5 - 



o.o L 



Xe 
Kr 
Ar 
I Ne 



s p d f g h i j 

Orbital symmetries 

FIG. 10: The second-order correlation energy when orbitals 
upto a particular symmetry are included in the virtual space. 



f 2 i 

collaborators [33] provide cumulative EcoIy for Ne up to 
i symmetry. Their converged result is -0.3836, this com- 
pares well with our result of -0.3811 calculated with or- 
bitals up to i symmetry. However, in our calculation, we 
get converged results of -0.3830 after including j and k 

(2) 

symmetry orbitals. For Ne, -Ecorr is considered converged 
when the change with additional symmetry is below mil- 
lihartree. However, for Ar, Kr and Xe orbitals of I and 

(2) 

higher symmetries are essential to obtain -Ecorr converged 

(2) 

up to millihartree. Since the magnitude of Ecorr increases 
with Z, a representative measure of convergence is the 
percentage change. The contribution from k symmetry 
to E { c llr are -0.0017, -0.0083 and -0.0176 for Ar, Kr and 
Xe respectively. These are larger than that of Ne, which 
is -0.0007. However, these correspond to 0.24%, 0.45% 



TABLE II: Cumulative second-order correlation energy when 
orbitals upto a particular symmetry are included in the virtual 
space. All the values are in atomic units. 





Symmetry 


Ne 


Ar 


Kr 


Xe 




s 


-0.0194 


-0.0210 


-0.0236 


-0.0247 




P 


-0.1920 


-0.2043 


-0.2479 


-0.2687 




d 


-0.3216 


-0.5401 


-0.9512 


-1.0419 




f 


-0.3589 


-0.6330 


-1.5213 


-2.2972 




9 


-0.3732 


-0.6695 


-1.7077 


-2.6879 




h 


-0.3786 


-0.6830 


-1.7843 


-2.8520 




i 


-0.3811 


-0.6891 


-1.8179 


-2.9238 




J 


-0.3823 


-0.6921 


-1.8343 


-2.9591 


k 


-0.3830 


-0.6938 


-1.8426 


-2.9767 



and 0.59% for Ar, Kr and Xe respectively, these com- 
pare very well with that of Ne 0.20%. For Ar there is 

(2) 

a variation in the previous values of -Ecorr, these range 
from the lowest value of Clementi [34] -0.790 to that of 
Ishikawa [32] -0.6981. Our value of -0.6938 is closer to 
that of Ishikawa. 

There is a pattern in the change of the correlation en- 
ergy with symmetry wise augmentation of the virtual or- 
bital set. There is an initial increase, reaches a maximum 
and then decreases. The maximum change occurs with 
the addition of p, d, d and / symmetry for Ne, Ar, Kr and 



Xe respectively. The pattern is evident in Fig 10 which 

(2) ■ 

plots the change in £corr with symmetry wise augmen- 
tation of the virtual space. The pattern arise from the 
distribution of the contributions from each of the core or- 



bitals. From Eq.( 19 1, depending on the core orbital com- 
bination ab, there are two types of correlation effects. 
These are inter and intra core shell correlations corre- 
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sponding to a — b and a b respectively. Among the var- 
ious combinations, the 2p 3 / 2 2p 3 / 2 , 3p 3 / 2 3p3/2, 3d 5 / 2 3d 5 / 2 
and 4d 5 / 2 4:d 5 / 2 core orbital pairings have leading contri- 
butions in Ne, Ar, Kr and Xe respectively. Here for Ne 
and Ar the leading pairs correspond to the valence shell 
but it is the last d shell for Kr and Xe. This correlates 
with the pattern observed in the symmetry wise augmen- 
tation. 



B. Third-Order Correlation Energy 

(3) 

We calculate the third order correlation energy Ecorr 
from the linearized CCSD equations. This is possi- 
ble when the first order MBPT wave operator fiW is 
chosen as the initial guess and iterate the linearized 
coupled-cluster equations once. The two-b ody wave op- 



(2) 

erator so obtained is fl 2 and from Eq.(21l E, 



(^o|^2^4 2; |^o)- This approach, however, is not appli- 
cable beyond third order. The reason is, starting from 
the fourth order correlation energy the triples contribute 
to AE and triples are not part of the linearized CCSD 

(3) 

equations. The results of Eeorr, obtained from our cal- 
culations, are listed in Table |TTT| For comparison, results 
from previous works are also listed. For Ne, Jankowski 
and Malinowski [39 reported a value of 0.0024. Their cal- 
culations were done with a limited basis set and hence, 
could leave out less significant contributions. The re- 
sults of Lindgren and collaborators [33] 0.0035 is perhaps 
more accurate and reliable on account of larger basis set. 
They include virtuals up to i symmetry and then ex- 
trapolate. Similarly, in our calculations we include vir- 
tual orbitals up to i symmetry, then extrapolate up to 

(2) 

k symmetry based on Ecott results. We obtain 0.0019, 
which is in better agreement with the result of Jankowski 
and Malinowski [39 . As expected, Ecorr increases with Z 
and to our knowledge, our results of Ar, Kr and Xe are 
the first reported calculations in literature. Interestingly, 

(3) 

Ecott is positive for Ne, Kr and Xe but it is negative for 
Ar. 



n(3) 



TABLE III: Third-order correlation energy in atomic units. 



Atom 




E :i 






This work 




Other work 


Ne 


0.0019 




0.0035" 








0.0024 6 


Ar 


-0.0127 






Kr 


0.0789 






Xe 


0.1526 







a Reference[33], 
i, Reference|39l. 



C. Coupled-Cluster correlation energy 

The MBPT correlation energies Econ converges with 

(2) 

relatively large basis set. For example, the Ecorr of Ne 
converge when virtual orbitals up to k symmetry are in- 
cluded in the calculations. This correspond to a total of 
224 virtual orbitals. Similar or larger number of virtual 

(2) 

orbitals are required to obtain converged Ecorr of Ar, Kr 
and Xe as well. However, it is not practical to have such 
large basis sets in relativistic coupled-cluster calculations. 
The n\ri? c , where n v and n c are the number of the vir- 
tual and core orbitals respectively, scaling of arithmetic 
operations in CCSD render computations with large n v 
beyond the scope of detailed studies. Hence, in the CCSD 
calculations, the size of the virtual orbital set is reduced 
to a manageable level and restrict up to the h symme- 
try. To choose the optimal set, after considering the most 
dominant ones, the virtual orbitals are augmented in lay- 
ers. Where one layer consists of one virtual orbital each 
from all the symmetries considered. 

The CCSD correlation energies for two basis sets are 
listed in Table ITVl The first is with a basis set consid- 
ered optimal and manageable size for CCSD calculations 
after a series of calculations. Then the next is with an 
additional layer of virtual orbitals. The change in the lin- 
earized CCSD AE with the additional layer of virtual or- 
bitals are 0.2%, 1.7%, 6.0% and 5.0% for Ne, Ar, Kr and 
Xe respectively. Changes of similar order are observed in 
the corresponding AE of the non-linear CCSD calcula- 
tions. It must be mentioned that, though the difference in 
AE is small, the computational cost of non-linear CCSD 
is much higher than the linearized CCSD calculations. 
The percentage changes indicate the basis size of Kr and 
Xe are not large enough. The orbital basis of Xe, with 
the additional layer, consists of 17 core and 129 virtual 
orbitals. This translates to ~ 5.0 x 10 6 cluster ampli- 
tudes, which follows from the n^n 2 c scaling of the number 
of cluster amplitudes. At this stage, the computational 
efforts and costs far out weight the gain in accuracy. To 
account for the correlation energy from the other virtual 
orbitals, not included in the CCSD calculations, we resort 
to the second order correlation energy. For this we calcu- 

(2) 

late Ecom- with the basis set chosen in CCSD calculations 

(2) 

and subtract from the converged Ecott- The estimated 
AE in Table ITVl is the sum of this difference and CCSD 
AE. This includes the correlation effects from i, j and k 
symmetries as well. For Ne, the estimated experimental 
value of correlation energy lies between the range 0.385 
and 0.390 [36j |40] . Our coupled-cluster result, estimated 
value, is in excellent agreement. 

The contributions to the correlation energy arising 
from the approximate triples are listed in Table [Vj As dis- 
cussed in Section |lIIB| the correlation energy diagrams 
corresponding to the approximate triples are grouped 
into three classes. Out of these we have selected a few: 
eight from 2p-2h and two each from 3p-lh and lp-3h. In 
Table |V] AE arising from these are listed. It is evident 
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TABLE IV: Correlation energy from coupled-cluster. All the 
values are in atomic units. 



Atom 


Active Orbitals 


AE( 


ccsd) 






Linear 


Nonlinear 


Ne 


17sl0pl0d9/9cf8/i 
18sllplldl0/10ff9/i 
Estimated 


-0.3783 
-0.3805 
-0.3905 


-0.3760 
-0.3782 
-0.3882 


Ar 


17sllplld9/9ff9/i 
18sl2pl2dl0/10pl0ft 
Estimated 


-0.6884 
-0.7001 
-0.7258 


-0.6829 
-0.6945 
-0.7202 


Kr 


22sl3plld9/9#9/i 
23sl4pl2dl0/10sl0ft 
Estimated 


-1.5700 
-1.6730 
-1.8480 


-1.5688 
-1.6716 
-1.8466 


Xe 


23sl4pl2dl0/105l0/i 
24sl5pl3dll/llflllft 
Estimated 


-2.5500 
-2.6874 
-2.9973 


-2.5509 
-2.6881 
-2.9979 



from the table, the contribution from lp-3h and 3p-lh 
are negative and adds to the magnitude of AE. Whereas, 
the contribution from 2p-2h is positive and reduces the 
magnitude of AE. 



TABLE V: Correlation energy arising from the approximate 
triples in the coupled-cluster theory. All the values are in 
atomic units. 



D subsumes contributions arising from Dirac-Fock 
and random phase approximation. We can thus expect 
this term to have the most dominant contribution. This 
is evident in Table fVT] which shows that the contribu- 
tion from T 1 - 1 ^ D is far larger than the others. The next 
leading term is DT^\ This is attributed to the 

larger values, compared to T^°\ of cluster ampli- 
tudes. The dipolc polarizability calculations with rela- 
tivistic coupled-cluster theory involve two sets of cluster 
amplitudes. These are the T(°) and cluster am- 
plitudes. As mentioned earlier, solving coupled-cluster 
equations is compute intensive. To minimize the compu- 
tational costs, we optimize the basis parameters «o and 
f3 to contract the size of the orbital basis set. 

One pattern discernible in the results is the better 

agreement between the T^ 1 ^ D results and experimen- 
tal data. The deviations from the experimental data are 
large when we consider the total (CCSD) result. For 
Ne, the deviation from experimental data is 2%, where 
as it is 9% in the case of Xe atom. We attribute this 



to the approximation in Eq. ( 50 ) and partly to the basis 



set. To confirm this, however, requires detailed computa- 
tions with higher order terms in the unitary transforma- 
tion. This poses considerable computational challenges 
and shall be addressed in future publications. We also 
implement the approximate triples excitation of the per- 
turbed cluster amplitudes and contribution to a are listed 
in the table. 



Atom 


Basis size 




AE 


















2p-2h 


lp-3h 


3p-lh 


TABLE VI: Dipole 


polarizability of the e 


round state of neu- 


Ne 


18sllplldl0/10 5 9/i 


0.00672 


-0.00145 


-0.00164 


tral rare-gas atoms (in a. u.). 








Ar 


18sl2pl2dl0/10pl0/i 


0.00805 


-0.00066 


-0.00192 


Contributions 


Ne 


Ar 


Kr 


Xe 


Kr 


22sl3plld9/9g9/i 


0.01546 


-0.00171 


-0.00305 


T (^ DT (o) 


2.7108 
0.0771 


11.3330 
0.0486 


17.2115 
0.0429 


27.7427 
-0.1495 


Xe 


19sl5pl0d9/552/i 


0.02011 


-0.00148 


-0.00260 


T (^ DT (o) 
T i^ DT {0) 


-0.0703 


-0.8264 


-1.2721 


-2.3286 












-0.0004 


-0.0001 


0.0002 


0.0027 












Total(CCSD) 


0.0053 


0.2490 


0.0439 


0.0786 












2.7225 


10.8041 


16.0264 


25.3459 




D. Dipole 


Polarizability 




Approx. triples 2.7281 
Exp. values" 2.670±0.005 


10.7360 
11.070(7) 


16.0115 
17.075 


25.2974 
27.815 



One constraint while using perturbed coupled-cluster 
theory to calculate dipole polarizability is the form of D. 
It is a unitary transformation of the dipole operator and 
expands to a non terminating series. For the present cal- 
culations we consider the leading terms in D. That 
is, we use the approximation 



T (1) D « T{ 



D + DT^ + DT, 



-,(0) 



T. 



(i)t 



DT, 



(o) 



DT 



(0) 



(50) 



The ground state dipole polarizabilities of Ne, Ar, Kr 
and Xe calculated with this approximation are given 
in Table VI Among the various terms, the first term 



a 18 . 



VI. CONCLUSION 

We have done a systematic study of the electron corre- 
lation energy of neutral inert gas atoms using relativistic 
MBPT and coupled-cluster theory. Our MBPT results 
are based on larger basis sets consisting of higher sym- 
metries than the previous works. Hence more reliable 
and accurate. Our study shows that in heavier atoms 
Kr and Xe, the inner core electrons in d shells domi- 
nates the electron correlation effects. This ought to be 
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considered in high precision properties calculations. For 
example, the EDM calculations of Xe arising from nu- 
clear Schiff moment. The dipole polarizability calculated 
with the perturbed coupled-cluster show systematic devi- 
ation from the experimental data. However, the leading 
term is in good agreement. The deviations might de- 
crease when higher order terms are incorporated in the 
calculations. From these results, it is evident that the 
basis set chosen is of good quality and appropriate for 
precision calculations. 
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